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We study the interaction effect in a three dimensional Dirac semimetal and find that two competing 
orders, charge-density-wave orders and nematic orders, can be induced to gap the Dirac points. 
Applying a magnetic field can further induce an instability towards forming these ordered phases. 

The charge density wave phase is similar as that of a Weyl semimetal while the nematic phase is 
unique for Dirac semimetals. Gapless zero modes are found in the vortex core formed by nematic 
order parameters, indicating the topological nature of nematic phases. The nematic phase can be 
observed experimentally using scanning tunnelling microscopy. 

PACS numbers: 71.55.Ak, 71.20.-b, 71.45.-d 


Dirac semimetals are materials whose bulk valence 
and conduction bands touch only at certain discrete 
momenta, around which the low energy physics is de¬ 
scribed by gapless Dirac fermions with linear energy dis¬ 
persion. The two-dimensional Dirac semimetal is realized 
in graphene and has been studied extensively. The three- 
dimensional Dirac semimetals were predicted to exist in 
NasBi and Cd 3 As 2 [l|-[3[ and these predictions were con¬ 
firmed in the recent angular resolved photon emission ex¬ 
periments ill- The three-dimensional Dirac semimetal 
contains multiple copies of Weyl fermions and thus, they 
can exhibit non-trivial topology. Different from Weyl 
semimetals, the gapless nature of Dirac semimetals re¬ 
quires the protection of crystalline symmetries. As a con¬ 
sequence, by breaking some of these symmetries, Dirac 
semimetals can be driven towards other exotic states such 
as Weyl semimetals |ilI3 and axionic insulators miiii. 

In this letter, we consider the mass generation of 
a three dimensional Dirac semimetal with two Dirac 
points protected by rotational symmetry, such as the 
one realized in NasBi. Three different complex mass 
terms will arise when interaction is included in the effec¬ 
tive Hamiltonian of a three-dimensional Dirac semimetal 
NasBi. One complex mass is generated by charge den¬ 
sity wave (CDW) that involves inter-Dirac-cone scatter¬ 
ing and breaks translational symmetry. The other two 
complex masses come from nematic orders that break 
three-fold rotational symmetry (Cs) by involving excita¬ 
tions with different spins but within a single Dirac point. 
Within the mean field approximation, we map the phase 
diagram and find that intra-Dirac-cone interaction can 
spontaneously break rotational symmetry and drive the 
system into topological nematic phases. Electron charge 
distribution is identified for nematic phases, which can 
be directly detected by scanning tunnelling microscope 
(STM). We further discuss localized states in topological 
defects as a consequence of topological nature of nematic 


phases. We would like to emphasize that since gapless 
Dirac cones are protected by rotational symmetry, a gap 
opening by breaking rotation symmetry can lower the 
energy of a Dirac semi-metal. Thus, the presence of ne¬ 
matic phases is generic in rotational-symmetry protected 
Dirac semimetals. 

Let us start by describing our model. The low energy 
physics of NasBi is well captured by the k • p type of 
Hamiltonian density iLo(k) around the T point (H 


Hoik) 


/M(k) Ak^ 0 0 \ 

Ak_ -M(k) 0 0 

0 0 M(k) -Ak- 

^0 0 —Ak-^ —M(k) j 


( 1 ) 


up to the second order in k, where M(k) = Mq — Mik^ — 
M 2 ikl-\-ky). The bases here are |p+, t)^ \P-A 
), where for a basis \ol,(j), a = s,p± is the orbital in¬ 
dex and (T is the spin index. The above bases 

are also denoted as |^), ||), | — ^), | — |) based on the 
total angular momentum of each state. Mq, Mi, M 2 
and A are material dependent parameters. The part 
of iLo(k) that is proportional to the identity is not im¬ 
portant and has been neglected. The energy disper¬ 
sion is E(k) = ±>/M^(k) + A‘^k-^k- and two gapless 

points are located at Ki = ^0, 0, (—1)^ y^Mq/Mi^ , with 
i G {1,2}. The low energy effective Hamiltonian around 
Ki and K 2 can be expanded from m, and it is given 
by Ho = in the second quantized lan¬ 

guage, where 


Alik) = iCkp,s,tHk,l,p,tHk,l,s,i,Ckp,p,l, 

C/c,2,s,t, C/c,2,p,t, C/c,2,s,t, Ck,2,p,i)^, 

Ho = Ak^ao (8) Ls - Akyao (8) r4 + mikz)a^ (8) r 5 ( 2 ) 

k = ikx,ky,kz) is the momentum relative to the Dirac 
points Ki, mikz) = —2y/MoMikz and cl ■ ^ ^ creates an 
electron with a orbital and spin a at Ki P k. We also 
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denote Ck^i^p^^a as Ck^i^p^a for brevity, a^r^a are Pauli 
matrices characterizing spin, orbital and valley degree of 
freedoms. P matrices are defined as Pi, 2,3 = cri^ 2,3 ^ 'ri, 
r 4 = (Jo (8) T 2 and Ps = ao ( 8 ) T 3 . It is easy to see that they 
obey Clifford algebra {ri,rj} = 

We note that Hq is the minimal model for Dirac 
semimetals with time reversal (TR) symmetry and inver¬ 
sion symmetry. To describe the effective Dirac behavior 
of electrons near Ki, we keep only the linear terms in k. 
It should be emphasized that including other higher or¬ 
der off-diagonal terms cannot open a gap at Ki and K 2 
since two degenerate states transform differently under 
three-fold rotational symmetry. 

The fermionic field operator T can be thought of as 
four copies of Weyl fermions, two with left-handed chi¬ 
ralities and the other two right-handed. Here, we focus 
on the case with charge conservation and thus, the mass 
terms can only be formed by interactions of two Weyl 
fermions with opposite chiralities and therefore, there 
are two possible mass terms. The first one comes from 
two Weyl fermions with opposite chiralities at different 
momenta and K 2 ). This term breaks translational 
symmetry and corresponds to CDWs. Such a term can 
also be found in Weyl semimetals and is responsible for 
axion insulator phases [uHii- Since Dirac semimetals 
can be viewed as two copies of Weyl semimetals that 
conserve TR symmetry, the gapped phase due to CDWs 
should also be thought of as two copies of axion insulator 
phases which are related to each other by TR symme¬ 
try. The second mass term couples two Weyl fermions 
at the same momentum {Ki or K 2 ). Since the gapless 
nature of Dirac semimetals at a fixed momentum is pro¬ 
tected by Cs symmetry, it is natural to expect this mass 
term to break rotation symmetry but preserves transla¬ 
tional symmetry. This corresponds to a nematic phase. 
Therefore, these mass terms should be generated by the 
following order parameters: 

CDAV . Da,^,(7 — ^ ^/c,l,Q;,cr^^?2,/3,cr ^5 

nematic : Na,/3,Ki =< > • (3) 

On the other hand, possible mass terms should then anti¬ 
commute with Hq and there are only six of such terms: 
(ao( 8 )ri, (ao( 8 )r 2 , ai ( 8 ) Ts, 0^2 C) Ts, (a3(8)r2. 

Based on the above analysis, we identify all possible mass 
terms and introduce 

Kp,2 = Ks,2 = Ai - A 2 , 

= Ks,i = = A 3 , (4) 

where A^’s are generally complex: Aj = {j G 

1,2,3). 

To dynamically generate these mass terms, we consider 
an effective interaction between different species of Dirac 


(a) (b) 




FIG. 1. (a) Phase diagram of interacting 3D Dirac semimetal 
NasBi. (b) In the nematic phase, the ratio between A 2 and 
Ai is plotted as a function of V/U. 


fermions as given by 


Hint = U EE Pi(k)pi{k) + V EE Pi{k)pj{k),{5) 

k i k i^j 

where pi = are the density operators. 

Here, the U term describes the interaction between two 
electrons near one momentum Ki while V term describes 
that of electrons between Ki and K 2 . This effective in¬ 
teraction can be obtained from the Coulomb interaction, 
as shown in the Supplementary Materials [14 1. 

The full Hamiltonian can then be treated within the 
mean field approximation (see the Supplementary Mate¬ 
rials [ 1 ^ for details), the free energy at zero temperature 
is then given by 


F=K 


MF 


E 

E^fcGoccupied 


.E,(|Ai|,|A2|,|A3|,0), (6) 


Here is the excitation spectrum in the mean field level, 
whose detailed expression is shown in the Supplementary 
Materials [l^ . 0 = O 1 —O 2 represents the phase difference 
between the two nematic order parameters. To minimize 
the free energy, a state where 0 = ^ is favored. We estab¬ 
lish self-consistency equations to map the phase diagram 
in Fig. dja). The semimetallic phase is relatively stable 
under weak interaction because the density of states van¬ 
ishes at Dirac points. As the interaction strength exceeds 
critical value Uc (W), tbe system develops a gap. In the 
large U (V) limit, the system favors nematic (CDW) or¬ 
dering. Starting from the bi-critical point {Uc^Vc), the 
system will go across a first-order phase transition at the 
U = V line between the CDW and nematic phases. 

The ordered phase of CDW is similar to that in Weyl 
semimetals, the physical consequence of which has been 
discussed in details in [IMl. For our system, the 
CDW is along the kz direction with the wave vector 
Q = 2^Mo/Miz, which can be in principle observed 
in an STM. Chiral modes have been proposed to exist 
at the dislocation line of Weyl semimetals, but since our 
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FIG. 2. The energy dispersion from realistic k • p theory and 
LDOS for one Bi layer, (a) and (c) are for free semimetal, 
while (b) and (d) are for the interacting case. In the LDOS 
plot, red (blue) color represent a large (small) LDOS. 


TR invariant system is a copy of two Weyl semimetals, 
we have two copies of chiral modes that are TR partners 
and thus, our system exhibits helical modes. 

What is really unique in the Dirac semimetals is the 
nematic phase. This nematic phase is actually a super¬ 
position of two inequivalent nematic orders Ai and A 2 
with a phase difference of By applying TR operation 
0 = 0^1 G ia 2 G) To, we find Ai breaks TR symmetry 
while A 2 preserves TR symmetry. In Fig. [I](b), the co¬ 
existence of two nematic orders is numerically confirmed. 
As the ratio V/U increases from 0 to 1 , we find that the 
ratio A 2 /A 1 decreases from 1 to 0. This indicates that 
the system spontaneously breaks TR symmetry in the 
nematic phase. Next, we will discuss several physical 
phenomena of nematic phases, which can be observed in 
experiments. 

The first observable is the charge distribution. Since 
the mass term of nematic orders couples |±|) to |±^), we 
expect the charge distribution in one unit-cell breaking 
three fold rotation. Since the charge distribution cannot 
be extracted from the effective Hamiltonian, we consider 
a more realistic k • p Hamiltonian based on the first prin¬ 
ciples calculations. The method has been successfully 
applied to the construction of the effective Hamiltonian 
of topological insulator materials 0 , and we only de¬ 
scribe our procedure briefly here. The eigen wave func¬ 
tions 8it k = 0 can be expanded in term of plane waves 
in the first principles calculations. Here, 36 bands are 
taken into account, denoted as |n) (n = 1,2,-•• ,36). 
We act the crystal Hamiltonian with periodic potential 
on the basis and obtain the k • p Hamiltonian = 

■ Pnm, where En is the eigen- 
energy for the n band at /c = 0 , m is electron mass 


and Pnm = {n\p\m) is the matrix element. We diago¬ 
nalize this 36 X 36 Hamiltonian and the energy disper¬ 
sion is shown in Fig. which qualitatively fits to 

that from the first principles calculations. In particular, 
a level crossing between conduction and valence bands, 
which corresponds to Dirac points, can be seen along 
the V — Z line. From the eigen wave functions, one can 
show that the conduction and valence bands indeed be¬ 
long to the I ± |) and | ± |) states, respectively. Thus, 
these two states cannot be coupled to each other along 
the V — Z line. As discussed above, the interaction can 
introduce the coupling between these two states in the 
nematic phase. Therefore, we can add a constant cou¬ 
pling between the | ± and | T f) states near the Fermi 
energy in our k • p Hamiltonian, leading to a gap open¬ 
ing, as shown in Fig. Efb). To show that the obtained 
states possess nematic orders, we calculate the local den¬ 
sity of states (LDOS) in the x-y plane for the Bi layer. As 
shown in Fig. Efc), without interaction, the maxima of 
the LDOS (red color) appear as an isotropic ring around 
the position of Bi atoms, corresponding to the p± orbitals 
of Bi atoms. The spatial distribution of LDOS respects 
three-fold rotation symmetry. After adding the coupling 
term between ±|^) and T|f) states, the isotropic ring 
evolves into two peaks pointing a certain direction, thus 
breaking rotation (see Fig. Efd)). This corresponds 
exactly to the nematic phase. Such electron density dis¬ 
tribution can be directly measured through STM. 

The second phenomenon is the appearance of gapless 
modes in topological defects of the nematic phase, re¬ 
vealing the topological nature of this phase. Complex 
mass terms A = |A|e*^ in a Dirac system are highly non¬ 
trivial in the sense that their phases 0 are identified as 
dynamical axion fields, which will give rise to bulk ax- 
ionic terms in the form of F^^Fpo- H, 

In 2D Dirac systems, complex mass terms will show up as 


22 , 


a [/(I) or Zn vortex structure in both graphene _ _ _ 

and TT-fiux square lattice [ 3 , [HI in the presence of in¬ 
teractions. As a consequence, zero modes will localize 
at the vortex centers carrying fractionalized charges. In 
3D Weyl/Dirac systems, those zero modes extend to ID 
chiral modes that go through the center of the vortices 
along the 2 ;-direction [ll|, ll2| . These are known as axion 
strings. As is in the case of CDW, we expect a similar 
physics to occur in the vortex of nematic order parame¬ 
ters. 

By applying in-plane vortex structures for the complex 
nematic order parameters, our system at fixed kz can be 
directly mapped into previous 2D Dirac systems. There¬ 
fore, zero modes are expected to show up at both Ki and 
K 2 . To verify this, a numerical calculation is performed 
in a layered 2 D vortex configuration. We keep the peri¬ 
odicity in the z direction, while placing open boundary 
conditions in the x-y plane. For simplicity, on a 32 x 32 
square lattice, we place a U{1) vortex-antivortex pair 
structure instead of the actual Z 3 vortices. We visual- 
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FIG. 3. (a) The fermionic spectrum in a U{1) vortex- 

antivortex configuration on a 32 x 32 square lattice with open 
boundary conditions. We choose the following set of param¬ 
eters: Mo = -0.6, Ml = -0.3, M 2 = -0.4, A = 0.4, |f/| | Ai | = 
0.25, |f/||A 2 | = 0.1. Gapless energy bands in red (green) are 
localized at the vortex (anti-vortex) center, (b) LDOS at 
Ef = 0 is plotted which clearly shows zero modes are local¬ 
ized at vortex or antivortex center. The red (green) dot shows 
the location of a vortex (antivortex) center while yellow (red) 
color represent a large (small) LDOS. 


ize these vortex structures in Fig. [3jb) by the arrow 
indicating phase information of following site-dependent 
order parameters |24| : 


Ai{x,y,k^) = |Ai 
A2{x,y,h) = k 


{lO — UJi){lO — LO 2)* 
\{u: - u)i){u! - UJ2)\'' 
A 2 I 


Ai 


■Ai{x,y,k^). 


( 7 ) 


Here, uj = x-\-iy is a complex coordinate and Uj = Xj-\-iyj 
is the complex coordinate of vortex {j = 1) or anti-vortex 
center {j = 2). As shown in Fig. [3fa), the bulk disper¬ 
sion is gapped while gapless modes penetrate the bulk 
gap twice at two different momenta. In Fig. [3jb), we 
plot the LDOS at Ep = 0 together with vortex config¬ 
urations in real space. It is confirmed that these modes 
are highly localized at the vortex (anti-vortex) center. 
The gapless nature of these modes relies on the fact that 
they are separated at different momenta, and requires 
the translational symmetry along the z-direction. 

So far, we have discussed the effects of interaction in 
driving Dirac semimetals toward other phases. How¬ 
ever, those phases along with their novel physical phe¬ 
nomena can only be realized under relatively strong in¬ 
teraction. To overcome this difficulty, one can apply a 
magnetic field along the z direction such that Landau 
levels emerge. Similar strategies have been applied to 
achieve quantum Hall ferromagnetism in graphene sys¬ 
tems [26|, I 27 I, where spin orbital coupling (SOC) is al¬ 


most absent. The strong SOC in Dirac semimetals, how¬ 
ever, tends to tilt spins. As a result, CDW and ne¬ 
matic phases are more likely to be favored than fer¬ 
romagnetism in Dirac semimetals. The Landau levels 
in Dirac semi-metals have been observed experimentally 



FIG. 4. (a) Landau level dispersion along kz with magnetic 
field B = lOT. (b) Phase diagram of NasBi under magnetic 
field and interaction. 


[isl-l^. Even though the higher Landau levels of NasBi 
are gapped, the lowest Landau levels (LLLs) are gapless 
at Ki {i = 1 , 2 ), see Fig. 01 We identify this degen¬ 
eracy to be a crossing between |5,t) and |p, |) states, 
which is protected from developing a gap by C 3 symme¬ 
try along (001) axis. To describe the low energy physics 
of the gapless LLLs, we define a four-component spinor, 

= (4,i,«,t’4,i,p4’4,2,s,t>4,2,p4)- 

(|4]) are reduced to: (1) Density Wave: Di = D 2 = 

(2) Nematic: Ni = N 2 = Ns^p^K^- 

Through a similar mean field analysis (see the Supple¬ 
mentary Materials), the free energy at zero temperature 
is given by F = Hmf - Ei Efe, \A«(FF+T where 

m{kz) = —2yJMi{Mo — ^^kz- ^ 1,2 are functions of or- 
der parameters Di ^2 and A^i, 2 , whose detailed expressions 
are explicitly shown in the Supplementary Materials Hi- 
By minimizing the free energy, we obtain the phase 
diagram as shown in Fig. Hb). Instability happens for 
arbitrarily weak repulsive interaction [l 2 |, ll3| and as one 
tunes the interaction to go across V/U = 1, the system 
undergoes a phase transition from a CDW phase to a 
nematic phase or vice versa. Let us focus on the nematic 
regime (Di ^2 = 0 ) and the corresponding self-consistent 
equations can be solved analytically. As is shown in the 
Supplementary Materials [ 1 ^, the critical temperature 
that characterizes a phase transition from semi-metallic 
phase to the nematic phase is 


^ 2e^VfK h 

Tc = - U eBS 

irks 


(8) 


where A is the momentum cut-off and Vf = \ | is 

the Fermi velocity. 7 = 0.577... is the Euler constant 
and ks IS the Boltzmann constant. We have consid¬ 
ered a sample with a finite area S in the x-y plane. 
When T < Tc^ non-zero nematic ordering will always be 
formed for arbitrary U. In the zero temperature limit, 
the magnitude of order parameter can be solved Umi: 

|A^i| « 7^6 This expression indicates that 

a larger energy gap will show up for a larger magnetic 
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field. This instability under magnetic fields is a direct re¬ 
sult of the finite Landau level degeneracy. This suggests 
the necessary condition for the instability is that the cy¬ 
clotron length is much smaller than the sample size. In 
the Supplementary Materials Hi , we further discuss the 
existing experiments studying LLs of Dirac semimetals, 
and predict possible evidence of nematic phases in STM 
measurements of Landau levels. 
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Supplementary Materials for “Nematic phase of Dirac semimetal” 


Microscopic derivation of Hint from Coulomb interaction 


In this section, we give a microscopic derivation of the interacting term Hint from the well-known Coulomb inter¬ 
action, 


Hcoulomh — 

k,k' ,q a,/3,cr,(T' 


(9) 


Here, V{q) = 27Te^/\q\‘^. a and (3 are orbital indices, while a and a' are spin indices. Then we apply a mean field 
treatment to Hcoulomb, 


Hcoulomb — 

k,k' ,q Q;,/3,cr,cr' 


■“ ^ ^ ^i^)[^k-\-q,a,a^k',/3,a' (c^+g^Q,^cr^^',/5,cr') + 

k,k' ,q Q;,/3,cr,cr' 

^[^k'-q,l3,a'^k,a,a ~ + i^k'-q,l3,a'^k,a,a)] 

^ ^ ^ ^(^)[{^k-\-q,a,a^k',P,a') X (c^/_g^/3^cr'^/c,a,cr) 

k,k' ,q a,(3,a,a' 


^^k' —q,^ ,a' ^k-\-q^OL,(T^^' ->3 1 ^' ^^k-\-q,cx.,a^^' ^*^'')^k' —q,l3 ,a' ^k,Oi,a \ • 


( 10 ) 


Naively, we are particularly interested in the scattering process between Weyl fermions with opposite chirality. As 
shown in the main text, we have identified all possible mass terms (order parameters): 


CDW : < c 
Nematic : < c 


\:-\-q,a,a^^^33,cr >— X Sk^k' —q—2Ko, 

l+q,a,C,0,l >= ^a,l3,Ki X 5k,k'-q- 


( 11 ) 


In the definition of nematic order N^^is^Ki, we have defined that both k and k' are effective momenta relative to bulk 
Dirac point i^ 2 =i, 2 - Based on Eq. (HID, we are ready to decompose Hcouiomb into different channels Hcouiomb = 
HcDW T H]\[Qjnatic' 


HcDW = E E 1^(3)(I-^«,/3,<t|^ X 6k,k'-q-2Ko “ -Da,/3,<T X 5k,k'-q-2Kocl^q,a,a^k’,/3,a 

k,k' ,q Q;,/3,cr 

Da,^,a X ^/c,/c'—g—2Ko^I'—g,/3,cr^^5Ck:5^) 

~ y{Q){\^a,P,a\ — ^a,/3,cr X cl._^^^^^^Ck-\-q-\-2Ko,/3,(T ~ X cl_^2Ko,f3,a^k,a,a) 

k,q a,(3,a 

~ y{Q){\^a,P,a\ — ^a,/3,cr X c\.^^^^Ck-\-2Ko,^,o- “ X C^k^2Ko,^,cr^k,a,cr) 

k,q a,(3,a 

~ (|^a,/3,cr| “ ^a,l3,a ^ ^\^a,a^k-\-2Ko,/3,a ~ Da,l3,a X C^k+ 2 Ko,l 3 ,a^k,a,a) ‘ (12) 

k a,13,a 


Here, k is the effective crystal momenta relative to the Dirac point (0, 0, — i^o)^ therefore k H- 2Kq is in the vicinity of 
the other Dirac point (0, 0, i^o)- When discussing the nematic phase below, k is the effective crystal momenta relative 
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to the Dirac point = (0,0, (—depending on the i index of fermionic operator Ck^i^a^a- 

-f^Nematic ~ ^ ^ ^ ^ ^ ^ [|-^Q,/3,i | ^ ^k,k'—q ^a,/3,i ^ 

i=l,2 /c,/c' ,q a,(3 

= EEE 

i=l,2 k,q a,/3 

= EEE 

i=l,2 k,q a,/3 

~ ^ ~ ^a,l3,i ^ ^k,i,a,t^k,i,P,i ~ X . 

2=1,2 k a,l3 


(13) 


In the above expressions, we have defined a CDW (Nematic) interaction strength U (V). It is interesting to notice 
that U and V take the same value 27re^/|g'p. In our phase diagram of mean field theory, U = V corresponds to 
the critical line separating CDW phase with nematic phase. However, in realistic materials, we expect one of the two 
phases will be favored, depending on the material details, which is beyond the scope of this paper. 


MEAN FIELD THEORY 

Starting from Coulomb interaction, we have shown that the essential physics is captured by inter-Dirac-cone scatter¬ 
ing (Hcdw) and intra-valley-scattering (i^Nematic)- This inspires us to write the effective density-density interaction 
Eq. (5) in the main text: 

Hint = U EE Pi(k)pi{k) + V EE Pi(k)pj(k), (14) 

k i k i^j 

where pi{k) = a^ki a a^k,i,a,a are the density operators. This effective interaction term is equivalent to both Eq. 
([12]) and Eq. ([T3|) . while illustrating the physics in a better way. Based on the form of order parameters, we could 
put constraints to the indices and further simplify the density-density interaction to be 

k a,P,aj^a' 

PiPj — ^ ^ ^k,i,a,a^k,i,a,ac\ji^^^CkJ,p,a‘ (15) 

k a,l3,cr 

Applying a similar mean field analysis to our earlier discussion, the interaction terms can then be written as 
Pipi = '^'^{\Na,p,Ki\‘^ - iVa,/3,Ki4,i,/3,4,Cfe4.a.t “ 

k a,13 

PiPj ~ (|T^a,/3,cr| — Da,p,crc\ 2^l3,cr^k,l,a,cr ~ T)a,/3,cr^I,l,a,cr^^52,/3,cr), (16) 

k a,P,a 

where the order parameters are defined in Eq. (3) of the main article. 

Then, the mean field Hamiltonian is readily obtained 

k 

„ ^ pii 

H22)' 

Hmf = Y1 E U{\N„,ppf + \N„,p,2f) + V{\D„,pM^ + \D^,p,^\^), (17) 

k a,i3=s,p 


where 


T(^) - (c/c,l,s,t5 ^/c,1,S,45 ^/c, 1,P,45 ^/c,2,S,t5 ^A:,2,p,t5 C/c, 2,S,45 ^/c,2,p,4) 5 


(18) 
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and Hint is an 8 x 8 matrix with each Hij to be a 4 x 4 block; 


Hit = 


Hi2 = V 



m(k) 

Ak^ 

-17Ar*,,i -UNlpA 


Ak- 

—m{k) 



-UN,^s,i 

-UNy^S,l 

m{k) —Ak- 

v- 

-UNs,p,i 

UNp^p^i 

1 

+ 

1 



— D* 
^s,p,-r 

o 

o 


V 

-d;,sa 

— D* 

0 0 


0 

0 

_ n* — n* 



1 0 

0 

_ n* — n* / 



(19) 


Since we are only interested in mass terms that can gap the system, we would like to only keep mean field terms 
that anti-commute with the original Hamiltonian: 


n: 


N*. = 0 


D:,P,a = = 0 

= Ks,l = Ai + A 2 

Kvfi = Ks,2 = Ai - A2 
D:,s,t = = As 

D;,pA = Kp,i = -^3 ( 20 ) 

Here, Ai (A 2 ) is the nematic order that spontaneously breaks (preserves) TR symmetry and breaks three-fold rota¬ 
tional symmetry. A 3 is the charge density wave order parameters that breaks translational symmetry. Also notice 
that these order parameters are generally complex: Aj = |Aj|e*^J {j G 1,2,3). Then, we can write down Hint in n 
compact form: 

Hint ~ Hq H\^ 

Ho = Ak^ao (8) Ts - Akyao (8) r4 -h m{k)ao (8) Ts, 

Hi = [/|Ai|(cos6>i(ao (8) Ti - sinOiao (8) r2) + t/|A2|(cos6>2<a3 (8) Ti - sin6>2<a3 ^8) r2) 

+H|A 3 |(cos 6 > 3 (ai (8) Ts - sm0sa2 (8) Ts). 

The full Hamiltonian is then given by 

H = y"^\Ho - Hi)^ + Hmf, 


( 21 ) 


E' 

k 

HA. 


Hmf = 4( —)3[f/(|Ai|^ + |A2n + yiAsll- 

TT 


( 22 ) 


and is the volume of the sample and A is the momentum cut-off. The first term can be diagonalized analytically 
to yield the eigen-energy 

Ek = ±[C/2(|Ai|2 + |A2|2) + 1/2|A3|2 + AH+k. + m(fc,)2 

±2U\A2\VV^\A3\^ + U^\Ai\^cos‘^{0i-02)]i. (23) 

The above expression is the excitation spectrum that shows up in the free energy in the main text. 


ANALYTICAL PROPERTIES OF FREE ENERGY IN EQ. (8) OF THE MAIN ARTICLE 

Let us first show that why 0 = ^ is favored. Let us define 

e(k) = [/2(|Ai |2 -h |A 2 | 2 ) + Asp -h AH+k- + m{kf, 

f{k,0) = 2C/|A2|\/Y2|^3|2 + C/2|^^|2cos2 0, (24) 

such that the free energy can be written as 

F = ffMF- 2 ^J(k, 0 ), 
k 

J(k,0) = ye(k)+/(k,0) + ye(k)-/(k,0). 


(25) 
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Notice that Hmf is independent of 6>, and 


W 


1 

2 


1 

V'e(k) + /(k, 0) 


\/e(k) - /(k,6») 


] < 0 , 


(26) 


/(k, = /(k, ^ ~ f) ~ 2 t/y|A 2 A 3 |. So /(k, ^ will maximize J and thus minimize free energy F. So this 

condition constrains 0 = ^. 

Now we are ready to write down the self-consistency equations: 

A 1 1 /* .3. dJ 

~ 4U (2A)3 y ^5Ai’ 

^ 1 1 r 3, dJ 

' “ 4f/ (2A)3 J 

^ iy 

Here A is the momentum cutoff in the integration. The self-consistency equations can be solved numerically and the 
solution gives rise to the phase diagram in Fig. 1 of the main article. Analytically, they can also give us some hints 
on the shape of the phase boundary. After some manipulations, the first and the third equations in Eq. (|27|) are: 


By setting A^ 


1 

U 

1 


1 1 

4 (2A)3 

1 1 

4 (2A)3 


/ 

/ 




Ve + / 

y+y|t 

Ve + / 


1 

7 ^’ 


+ 


Ve- / 


0, we arrive at the critical interaction strength 


1 


1 


1 1 
2 (2A)3 


/ 




+ m{ky 


(28) 

(29) 


LANDAU LEVEL OF NA 3 BI AND THE SELF-CONSISTENT EQUATIONS 


Under a magnetic field that is oriented along the ^-direction, minimal coupling requires n = k ^A. Defining the 
magnetic length to be / = the commutation relation of tt is then given by 

, , ieB i 


'TTyJ — 


h 




(30) 


where we have chosen the gauge A = (0, Bx, 0). We can then define creation and annihilation operators in terms of 
TT as follows 


a = —=iT-, = —=71-1-, [a, a’^'l = 1. 

a/2 a/2 ^ ^ 


From the commutation relation, we find that 


2 2 , 1 

+ 7 r; = -p{aU+-). 


(31) 


(32) 


Then, by choosing the following trial wave-function \l/ = {fi 4 ’n, f 2 4’n-i^ fi4’NY', the Hamiltonian density 

can be written down as 




( M Afj^ 0 0 ^ 

Att- —M 0 0 

0 0 M —Att- 

^ 0 0 —Att^ —M j 

M+ fvW 0 0 \ 

fVM 0 0 

0 0 M+_1 ’ 

0 0 -j-V2N / 


( 33 ) 
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where 


M = Mo - Mikl - + i), 

Ci(pN = V^^AT-l* 

The Lowest Landau levels (LLL) are then given by = 0, 




/M(+ 0 0 0 \ 

0 0 0 0 
0 0 0 0 
VO 0 0 Mq J 


where 


M^k.) 


±MoTMiklT^- 


(34) 


(35) 


(36) 


This indicates that only the LLLs from |^) and | — f) states are gapless, with the gapless nodes located at Ki = 

(0,0,(-l)y^(Mo-^)). 

When considering instability problem of NasBi under strong magnetic field, the gapless LLLs are composed of the 
following states: 


\l) = \sA), \-l) = \p,i), (37) 

where we could define the following order parameters, 

Nematic : Ni = Ns,p,Ki, N 2 = Na,p,K 2 

Density Wave : Di = Ds^s,^:, -D 2 = -Dp.p.t (38) 

The Hamiltonian is then given by 

^ = E ^'^(-ffo +-ffint)^'+-ffMF, (39) 

k, 

where 

Ho = m(kz)Tz ® (Tj, 

/ 0 -UNI -VDl 0 \ 

-UNi 0 0 -VD*2 

- _yjj^ Q Q , 

\ 0 -VD 2 -UN 2 0 / 

Hmf = ^[t7(|iVi|2 + |7V2|2) + I/(|Di|2 + ID 2 I")]. (40) 

k 


Here, m(kz) = — ^)kz. The matrix part Hq + Hint can be diagonalized analytically, and the eigen- 

energy for occupied bands are —^m{kzY + 

6 = \[UWNi\^ + \N2f) + V\\D^\^ + ID 2 I") + y[/4(|7ViP - |iV2|2)2 + F4 (|Di|2 - |D2|2)2 

+27/21/2(1 ATiP + A2|2 )(|Di| 2 + |D2|2) + 8C/2y21ATiAr2DiD21 cos),).^, - 4>D2 - (t>N, + M 
6 = \[U^{\Nl\^ + \N2f) + + ID 2 I") - v/[/4(|/ViP - |iV2|2)2 + 1/4 (|Di|2 - |D2|2)2 


+2C/2l/2(|iVi|2 + |/V2|2)(|Di|2 + |D2|2) + 87/2^2 |iViiV2DiD2 | COs[,).Z,, - ct>D2 " 0iVx + (41) 

Since we are especially interested in the magnetic instability in the nematic regime, we can set density order 
parameters \Di\ = \D 2 \ = 0. Then in the mean field level, single particle Hamiltonian Hq + Hint have four energy 
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eigenvalues: = zband = ±Quantum partition function at finite 

temperature ksT = ^ {ks is the Boltzmann constant) is given by 

Z = Tre-^^ 

= ^(|iv,P+|iv.|^) ^ ^ e-^^f){l + e-^^i“)(l + e-^^^)(l + 

kz 

= g-/ 3 ^ti(|iV.|=+|iV.|=) X ^( 2 cosh^) 2 ( 2 cosh^) 2 . (42) 

kz 

Free energy E of this system is given by 

c = -iiogZ 

^ —(/(|A?iP + |A^2|^) - ■|^llog(2cosh^^) + log(2cosli^^)|. (43) 

kz 

Minimizing F with respect to |V| (i=l)2), we obtain the following self-consistent equations; 




2 IA 

TT 


ciw,i - 1 : tanh 

kz 


PEl U‘^\Ni\ 

~w~ 


Notice that the self-consistent equation for each order parameter is decoupled from each other. Therefore, we will 
discuss only one of the two nematic orders, for example Ni. 


FINITE TEMPERATURE EFFECT 


In this section, we will be discussing how a finite temperature will affect the appearance of different phases. In 
general, there should exist a critical temperature Tc that characterizes a finite temperature phase transition from a 
nematic (or CDW) ordered phase to an unordered gapless phase. At critical temperature Tc, order parameter vanishes 
so that we can perform the integration in the self-consistent equations: 


1 = 


U p dk^ 

aU-a“ 

/: 


tanh 


4A 


dkz tanh 
-A 

Ug 


I3E+ 1 
I3\m{kz)\ 


I 


2vfA 

Ug_ f 
VfA Jo 


dm{k 

l3vfAl2 


) tanh 


\m{kz)\ 

^\m{kz)\ 


i{kz)\ 


dx 


tanh: 


VfA 

Ug_ 

VfA 

Ug 


[(tanh X log x) — f 

Jo 


l3vfAl2 


dx- 


logx j 

cosh^ X 


[log 


PvfA 


f 


dx log X 


logx 
cosh^ 1 


[log - log J 

VfA 2 TT 

= Ea log 

VfA ^ ZirksTc ’ 

(44) 

where 7 = 0.577... is the Euler constant and N{0) is the density of states in ID. In the integration measure, we have 
considered the Landau level degeneracy in the x-y plane: 

_ _ eBS 

^ 27r/^ 


h 


( 45 ) 
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(a) (b) 




FIG. 5. Scaling of order parameter U\A\ (energy gap) are shown in: 
adopted the parameters from Ref. [2] and obtained Vf ~ 1.9 eV-A. 
calculations. 


(a) Area S = l/xm^ and (b) B = IT. Here, we have 
A momentum cut-off A = 0.2 is applied for the 


Here, S is the surface area of a Dirac semimetal sample spanned in the x-y plane. We also take the low temperature 
limit T ^0, so that oo. Therefore, we arrive at the relation between critical temperature Tc and interaction 

strength 1 /, 


^ 2e^VfA 2e^VfA h 

Tc= - ug = - u ess , 

irkB ttke 


(46) 


In this expression, we could clearly see that a larger U will naturally lead to a higher Interaction strength, however, 
is usually determined by the intrinsic properties of a material, and can barely be changed. Instead, we can increase 
the magnitude of the applied magnetic field which will enhance the transition temperature in a similar way. A simple 
estimation can be made for Tc'. if we take the sample in-plane area S = magnetic field 5 = 1 T, interaction 

strength U = 0.001 eV, then Tc turns out to be 1000 K. However, if sample area S is decreased to 0.5/im^, Tc = 210 
K. If sample area S is further decreased to 0.2/rm^, Tc = 1.8 K. Decreasing sample area is equivalent to decreasing 
magnetic field, since both quantities will influence Landau level degeneracy in the same way. This strong scaling 
behavior reflects the essential role of Landau level degeneracy in our discussions. Therefore, to observe the ordered 
phase (either nematic phase or CDW phase) we proposed, it is very important to prepare a sample of good enough 
quality and apply strong enough magnetic field. 


ZERO TEMPERATURE LIMIT AND GAP SCALING 


Next, let us look at the zero temperature limit. The free energy can be simplified to 

F = —[/(|iVi|2 + |iV2|2) - Y.{Et+E+). 


(47) 


Since |Vi| and IV 2 I are decoupled in the self-consistency equations, for |Vi| the self-consistency equation is given by 

dE 


5|Vi| 

2LK 


= 0 


TT 




U‘^\Ni 


J_ - J_ / 

Ug ~ 4A J_, 


dkz 


mik^Y + U‘^\Ni\^ 

1 1 


■log[ 


2vfA . 


(48) 


/-A \/m(fc^)2 -KC/2|jVi|2 27ru/A‘“'’'[/|Vi 

Here, A is the momentum cut-off and we define the Fermi velocity as vj = \ |. Then, the interaction-induced 


energy gap is 


|Vi 


2VfA 2nvfA 

» = 


2vfA h 

e 


U eBS 


(49) 
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FIG. 6. Dispersion of zeroth Landau levels are plotted for: (a) B = IT, |A| = 0.0 eV. (c) B = IT, |A| = 0.1 eV. The 
corresponding DOS plots are shown in (b) and (d). All parameters are adopted from Ref. [ 2 ]. 


Therefore, for an arbitrarily small 1/, a non-zero order (gap) will be developed. 

Based on Eq. (|49]), we are able to check the scaling relation of order parameter (gap) in terms of magnetic field 
B and sample area S. Numerically, these scaling relations are shown in Fig. [5l (a) We keep area S = l/rm^ and 
change the magnetic field B. (2) We keep B = 1 T and change the area S. Since Landau level degeneracy g ^ S x B, 
increasing either B or S will both increase the magnitude of interaction induced gap | A |. 


DENSITY OF STATES (DOS) AND POSSIBLE EXPERIMENTAL DETECTION 


To study possible interaction effect in a rotational symmetry protected 3D Dirac semimetal, we have proposed in 
the main text to visualize charge distribution by performing local density of states (LDOS) measurement with an 
STM setup. The appearance of an anisotropic charge distribution is identified as a key feature of the nematic phase. 
From a different perspective, the development of nonzero ordering also results in a finite gap in the energy spectrum. 
The energy gap of a system, however, is always ready to be read directly from the DOS measurement near the Fermi 
level, with the help of an STM setup. Therefore, in this section, we will discuss in details about the DOS feature 
of a Dirac semimetal sample placed in a strong magnetic field, where magnetic catalysis will assist the formation of 
ordered states. 

First of all, we would like to point out that a DOS measurement (or equivalently gap measurement) is only a direct 
evidence of the formation of a gap (symmetry breaking). However, such DOS measurement cannot distinguish a 
nematic phase from a charge density wave (CDW). Therefore an additional LDOS measurement is always necessary 
to determine the patterns of ordering before any conclusion can be reached. 

To start, we first consider a simplified problem where only the lowest Landau levels (LLLs) are present. The 
Hamiltonian of two gapless LLLs is = (Mq — Mik^ — Here (Jx,y,z are Pauli matrices defined under 

the bases |^) = (|^), | — |))^- In the discussions below, we will focus on the case of nematic phase where translational 
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E 


(f) 



FIG. 7. Dispersion of Landau levels in Dirac semi-metals are plotted for: (a) B = 20T, | A| =0.0 eV. (d) B = 20T, | A| = 0.01 
eV. The corresponding DOS plots are shown in (b) and (e). In (c) and (f), we zoom in to the red block region to get a better 
view of the DOS around E = 0. 


symmetry is preserved. Then a complex nematic order parameter A can show up in the off-diagonal part of Hq 


LLL 


H 


LLL 


(M = 


Mo 


Mikl 

A* 


M2 


-(Mo - Mikl - #) 


(50) 


Generally, for a one dimensional Hamiltonian H{kz)^ the DOS p{E) at energy E can be expressed in terms of retarded 
Green function G^{E^ 

k)- 1 

E-H{k,)+ir]’ 
p{E,k,) = -hm{Ti[G^{E,k,)]}, 

TT 

P{E) = I ^p{E,k,) = -hm J ^Tr[G^{E,k,)]. (51) 

Here. 77 <C 1 is a small number to avoid singularity. We have calculated both p{E^kz) and p{E) for and 

obtained band dispersions as well as the corresponding DOS figure. The DOS p{E) has an arbitrary unit because its 
calculated value is determined by the value of 77 we are choosing, and therefore only the relative magnitude of DOS 
within the same DOS plot is physically meaningful. As shown in Fig. [ 6 ] (a) and (b), when |A| = 0 , the system is 
gapless and the DOS of the ID Dirac point {E = 0) is finite. Notice that in Fig. [ 6 ](b), the DOS is diverging (peaks 
of DOS) at two different E^ which corresponds to two band extreme around E = ±0.1. When | A| = O.leV is turned 
on, the system is gapped (Fig. [ 6 ] (c)) and the DOS within the energy gap is suppressed. A new band edge formed 
around the energy E = —0.003 and E = 0.009, leading to two additional DOS peaks. 

Next, we consider a more realistic model where higher Landau levels are present (See Eq. (|33]) for details). As 
shown in Fig. [7] (b) and (e), if we tune E continuously, a peak of DOS will show up when E coincides with the band 
extreme of a Landau level. If we focus only on the low energy DOS around E = 0 , as shown in the red block regions 
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in (b) and (e), the DOS plots in (c) and (f) capture the main features of earlier discussions in Fig. [6l Experimentally, 
Fig. [3 (f) will be a smoking-gun signature of interaction induced ordering in rotational symmetry protected Dirac 
semi-metals. 

Unfortunately, the existing experiments are not intended for finding the nematic phase, although all the necessary 
conditions should already exist. The most closely related experiment is the STM measurement of Cd 2 As 3 from 
Yazdani’s group [28|. They even include a discussion of Landau level spectrum for magnetic fields along different 
directions (Fig. 4d and e in [HI), showing that two zero Landau levels will cross each other for magnetic field 
along [001] direction and anti-cross each other for [ 112 ]-directional magnetic field. Our prediction is that even for 
[001]-directional magnetic field, one still finds an anti-crossing behavior due to interaction effect. However, the STM 
measurement in Yazdani’s experiment is implemented on the [ 112 ] surface, which breaks C 4 rotation by itself. This 
prevents the observation of the nematic phase. To search for nematic phases, an STM measurement along the [001] 
surface is required. 

Another related experiment is the quantum oscillation of magneto-transport measurement in Cd 2 As 3 [? ]. Landau 
level splitting is resolved by rotating magnetic fields in the quantum oscillation measurements. However, this exper¬ 
iment can only reach the Landau level N > 2. Thus, to observe our prediction, one needs to further lower electron 
density to reach the truly quantum limit with experimentally feasible magnetic fields. In Ref. [HI, the quantum limit 
is reached at around 43T. However, the magnetic field is applied along the [112] direction, which again breaks the C 4 
rotation symmetry. 


